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Abstract 

We report a method for the efficient evaluation of analytic infrared (IR) intensities within gen- 
eralized Kohn-Sham density functional theory using Gaussian orbitals and periodic boundary 
conditions. A discretized form of the Berry phase is used to evaluate a periodic dipole moment and 
its derivatives with respect to in-phase nuclear coordinate displacements. Benchmark calculations 
are presented for one-dimensional chains of water molecules and poly(paraphenylenevinylene). 
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I. INTRODUCTION 



Theoretical prediction of vibrational spectra is a valuable tool for the characterization of 
compounds and materials which are either not synthesized or difficult to conduct experiments 
on. It requires at least two types of data: positions of peaks and their intensities. In the 
vibrational spectroscopy of periodic systems, peak positions are associated with phonon fre- 
quencies. Calculation of phonon frequencies from analytical second derivatives have recently 
become available for Kohn-Sham (KS) density functional theory (DFT) and Hartree-Fock 
(HF) methods, with plane- waves (PW^ and Gaussian type orbitals (GTO).- 1 ^ On the 
other hand, peak intensities have received less attention, even though they are necessary for 
comparison of theoretical spectra with their experimental counterparts. Part of the prob- 
lem is that only in the 90's was their correct theoretical expression derived.- 1 ^ Since then, 
first principle studies with KS-DFT and PWs have become common for infrared (IR) and 
Raman spectra of solids.-^ However, we are aware of only one work on IR intensities with 
GTOs. Although GTOs have many advantages,- 1 ^ especially for orbital dependent func- 
tional, GTO-based techniques are mostly developed in the quantum chemistry community, 
and therefore, are underrepresented for periodic systems. Thus, in this work, we would 
like to present a simple scheme for evaluation of IR intensities for periodic systems within 
KS-DFT and HF with GTOs. 

The first order IR intensity of a fundamental transition exciting the jth normal mode is 
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where Na is Avogadro's number, c is the speed of light, d x _ z are the Cartesian components 
of the dipole moment, and Qj is the normal mode coordinate.— Since any normal mode is a 
linear combination of atomic Cartesian coordinates (Rj), IR intensities are calculated from 
dd/dHi. There are two main features arising for these derivatives under periodic boundary 
conditions (PBC). First, one should consider derivatives with respect to in-phase nuclear 
coordinate displacements 

oo 

dd/dR i = Y / dd/dR ig , (2) 

g=0 

where g is the index of an unit cell. Although in-phase vibrations do not represent the whole 



phonon spectrum of the periodic system, only these vibrations are IR active.— Second, under 
PBC, the dipole moment per unit cell cannot be straightforwardly expressed as a matrix 
element of the position operator jjj&iiLM According to the modern theory of polarization,-^ 
the dipole moment per unit cell is a geometric quantum phase or a Berry phase. 15 A non- 
trivial quantum phase usually appears in cases when the system under study is coupled to 
the rest of the universe through some external parameter The Berry phase 7 can then 
be written as 



where ip(C) is the wavefunction and C is a closed contour. In the case of periodic dipole 
moment, the parameter ( is represented by a reciprocal vector k running over the Brillouin 
zone. Therefore, within single-determinant methods, Eq. [3] necessitates evaluating d/dh 
derivatives of the Slater determinant. Since the latter is obtained in the self-consistent field 
(SCF) procedure which usually adds to the solution arbitrary complex phase prefactors 
depending on the k value, it is not trivial to differentiate the SCF solution with respect 
to k. Two main approaches to this problem have been proposed. The first method in- 
volves a discrete representation of <9/<9k derivatives through matrix elements at neighboring 
k points. This representation is invariant with respect to phase arbitrariness which could 
arise from the SCF procedure.- 1 ^ The second approach introduces phase factors which can- 
cel the SCF arbitrary phases and allow one to evaluate <9/<9k by employing regular finite 
differences techniques.— Here, we start from the periodic dipole expression arising in 
the first approach. We consider the first approach as a simpler and computationally more 
robust alternative, since it does not require introduction of complex phase functions and 
band resolution techniques. Therefore, it can be seen as a black-box technique that has 
only the number of k points in the Brillouin zone sampling as an input. In addition, this 
method has been tested with plane waves before,— >2i and its robustness was demonstrated 
not only for the linear response regime but for finite electric fields as well. Yet another 
approach to the periodic dipole derivatives is possible in the case of PWs, where these 
quantities require only dipole matrix elements between occupied and virtual Bloch func- 
tions (4>nk\ r \4>ak) ■ These "transition" dipole elements can be transformed into the following 




(3) 
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well-defined expressions^ 

(<M*\M ~ ea(k) _ en(k) > W 

where i^scF is the crystal SCF Hamiltonian, e a (k) and e n (k) are orbital energies of empty 
and occupied bands. However, in the case of GTOs this approach could be less beneficial 
since Pulay's type of terms arise. To the best of our knowledge, the only implementation 
of IR intensities with GTOs was reported by Jacquemin and coworkers 24 using the second 
approach. Thus, the focus of the current work is in adapting the discretized Berry phase 
approach for evaluation of IR intensities with GTOs. 

This work can be also seen as complementary to our previous work on zone centered 
phonon frequencies.- One should keep in mind that IR intensities are always obtained with 
vibrational frequencies, and therefore, any algorithm for IR intensities can be seen as an 
efficient one if its execution time is negligible with respect to that for vibrational frequencies. 



II. THEORY 

A. Berry phase approach to periodic dipole moment 

We first consider the formulation of periodic dipole moment through the Berry phase 
approach in localized basis sets proposed by Kudin and coworkers.— We restrict our con- 
sideration to closed-shell one-dimensional periodic systems aligned along the z axes. This 
makes our description simpler without losing essential details. Throughout this work we will 
be using Bloch orbitals 

1 N c M 

<Mr) = ^=EE ^)C, n {k)e iaka (5) 

V iV c 9= M =l 

expanded over Gaussian atomic orbitals (AO) u 9 (r). Here, small Greek letters are for AO 
indices with subscripts denoting the unit cell number. In Eq. (j3j), N c is the total number of 
unit cells in the system, C^ n (k) are crystal orbital (CO) coefficients from the SCF problem 
under PBC, 10 a is the length of the translational vector, and M is the total number of basis 
functions per unit cell. The cell-periodic part of a Bloch orbital corresponding to the nth 
band 

u nk (r) = e~ ikz (l) nk {r) (6) 
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can be expanded as 

N c M 
V iV c 9=0 M =l 

According to the Berry phase approach,-^ the z component of the dipole moment per unit 
cell is 

d =- l ™TlJ BZ dk J dvu* nh {v) d -^. (8) 

Here, we avoid introducing an additional z subscript, since unless otherwise stated we always 
will consider the longitudinal part of the dipole moment. Substituting u n k{v) in the dipole 
expression by Eq. ([7]) we obtain the following representation for the periodic dipole moment 

d = d l + d 2 , (9) 

occ. M „ 

i = EE dk c; n (k)c un (k) (io) 

x J2e ik °%Z%-gaS°l), 

9 



occ. M 



= - lm EE j BZ dk (n) 



9 



where 



fJ,V 



dr /j, (r)zu g (r), (12) 



and 



S°Jt = / dr to(r) Vg (r), (13) 



are dipole and overlap matrices in the AO basis. Calculation of the d\ term is more efficient 
to perform in the AO basis 

di = J2 P °A Z lt-9*S°t), (14) 

fJ.^,9 

where P®® is the electron density matrix. The matrix (Z^ — gaS^) can be further split into 
Hermitian 

ZH = 2% - f S% (15) 



and anti-Hermitian 



ga 



2 ^ 



(16) 



parts. Both matrices Z°^ u and Z®® are real, and we refer to Hermiticity of their reciprocal 
space counterparts Z^ v {k) and Z^ y {k). Since we contract both parts with the Hermitian 
matrix P^, the anti-Hermitian part can be safely omitted. This leads us to the following 
working expression 



(17) 



The d>2 term requires special care because CO coefficients usually contain arbitrary k- 
dependent phase factors from a diagonalization procedure. Therefore, one cannot apply 
regular numerical differentiation schemes directly to CO coefficients. To cope with this 
difficulty we use a discretization approach 2 - which starts from rewriting Eq. (fTTj) in the form 



d,2 = — Im / dk tr 

'BZ 



C\k)S(k)^ 



(18) 



where only the occupied orbital part of C{k) is engaged. For further discussion we introduce 
the following O-by-0 matrix 



E(k,k') = C\k)S(k)C(k') 



(19) 



where O is the number of occupied bands. Treating k' as a variable and k as a parameter 
we can write 



k=k' 



Keeping in mind the orthogonality relation E(fc, k) = 1 we rewrite Eq. [201 as 
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The trace operation commutes with the differentiation 
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(20) 
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Applying the well-known matrix relation tr \a{A) = In det(A) to the H(k,k') matrix we 
obtain 



_d 

dk 



7 tr [In E(fc,fc') 



d_ 

k=k> dk 
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7 ln [det E(Jfe, k') 



k=k' 



(23) 



Since In det E(fc, k) = 0, using simple rectangular discretization for the k' derivatives and 
the Brillouin zone integration we arrive to the following discretized form^ for d^- 

N k 

d 2 = -Im^ln det E(k j ,k j+1 ) (24) 
3=1 

= -ImJ^ln det C t (A; J )S(A; j )C(fc J+ i). (25) 

3=1 

By introducing arbitrary phase factors e ^ in front of CO coefficients, one can easily see 
that they cancel each other after the summation over the Brillouin zone. Please note that 
one cannot evaluate Eq. f|T8|) by using a perturbation (linear response) expression to build 
C(k + Ak) from C(k). As was pointed out in Ref. |y, a perturbation theory implicitly uses 
the so called "parallel-transport" gauge which necessitates phase equality between C{k) and 
C(k + Ak), and as a consequence produces zero Berry phase. 

Similar to the molecular case, the transverse components of the periodic dipole moment 
can be written as 



where q = x or y and 



Q% = J di7i (rV s (r). (27) 

Therefore, their differentiation is very similar to that of the d\ part of the longitudinal dipole 
moment, and it will not be considered further. 



B. Dipole derivatives 

Differentiation of Eq. ffTTl) with respect to nuclear coordinate displacements gives 

dR { ^ dR { ^ ^ dRi ' 1 ; 

Both terms in Eq. (128]) are calculated with standard techniques: the coupled perturbed 
SCF (CPSCF) procedure for the density derivatives,- 1 ^ and a regular integral derivatives 
evaluation technique. 27 The CPSCF procedure evaluates the same density derivatives as in 



the case of vibrational frequency calculations; therefore computational overhead of the d% 
part is negligible, since it only involves calculation of one-electron terms. 

In order to obtain the d<i derivatives we adapt the algebraic approach discussed in Ref. 
The main step can be summarized as application of the following algebraic identity 
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_d_ 
dR 



{In det[E(R)]} =tr 



jr'(R)A E(H ) 



(29) 



which is an analog of Eqs. [22] and [23j This was proven from a different point of view in 



Ref. 
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The differentiation of fcj+i) gives rise to two types of terms in the d.2 derivatives 



+ 4 2) 



dC 



dR, 1 \dRiJ 1 \dR, 
which correspond to Pulay's type of forces from overlap integral derivatives 

j=l m,n 
X ^ran(^) %+l) 
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and response of CO coefficients 
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4 2) = -Im ^^S m 1 „(A; j ,A; j+1 ) 
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(30) 



(31) 



(32) 



(33) 



dR, 



^fiv{ kj ) Cum {kj+i 



dRi 



To complete the construction of the dipole moment derivatives we will express CO derivatives 
via the response matrix U obtained in the CPSCF procedure^ 

dc; n (kj) 



2^ Una(kj)C* a (kj). 



Thus we arrive to 



Nk occ. vir. 

4 2) = -Im ^2^2^2^ n (kj, k j+1 ) 

j=l m,n a 



(34) 



(35) 



x \ U* a (kj)Y> am (kj , kj + i) 

~\~£'na{kji kj-\-\)U arn (kjj r i / 



TABLE I: Optimized geometry (A and degrees) of a one-dimensional water chain used to bench- 
mark the accuracy of our analytic method with different basis sets. See Fig. [I] for the definition of 
the geometrical parameters. 







HF 




PBE 


6-31G 


6-311++G** 


6-31G 


6-311++G** 


Roh 


0.9625 


0.9477 


1.0213 


0.9843 


roH 


0.9492 


0.9410 


0.9818 


0.9692 




1.9084 


2.1349 


1.6367 


1.9158 


a 


2.7459 


2.9547 


2.5709 


2.8048 


a 


108.58 


104.31 


105.77 


102.76 



Since we do not evaluate R« dependent phases, the application of the linear response Eq. ( 13411 
is acceptable. 

The computational complexity of the presented scheme is negligible with respect to that 
of vibrational frequency calculations. For the d 2 part we have several matrix multiplications 
of M-by-M matrices and one O-by-0 matrix inversion for each k point. 

Two- and three-dimensional generalizations of the presented method can be done in the 
same way as for plane-wave implementations; therefore, we refer the interested reader to 
Ref. 20 for details. 



III. BENCHMARK CALCULATIONS 
A. Algorithmic tests 

To validate and assess the present formalism, which has been implemented in the devel- 
opment version of the GAUSSIAN program,— we have tested it on a model one-dimensional 
chain of water molecules (see Fig. [1] and Tabled]) with the geometry optimized within the HF 
and PBE methods.— In order to describe hydrogen-bond interactions properly, we employ 
the 6-311++G** basis set in most of our tests.— 
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In Table |TT] we compare dipole derivatives obtained in periodic calculations with those 
from oligomeric estimations. Given that the 6-31 1++G** basis set is too large to perform 
vibrational frequency calculations for long oligomeric chains, we have used the 6-31G basis 
set for this comparison. The oligomeric results are calculated by the difference scheme which 
is less prone to edge effects 



where 




here, the indices i and m enumerate atoms within a molecule and molecules within an 
oligomer, respectively. Numerical periodic results are obtained by calculating finite differ- 
ences between unit cell dipole moments in different geometrical configurations. As expected, 
oligomeric calculations need more molecules to converge in the periodic direction than in 
the others. The same trend can be seen in analytic periodic calculations with respect to 
the number of k points. The numerical PBC calculation is less sensitive to the number of k 
points, however, they have a fixed error due to the finite difference scheme of differentiation. 
This suggests that the periodic dipole moment converges faster than its derivatives with the 
number of k points. 

As a matter of practical interest we illustrate the dependence of the five normal mode 
frequencies and their IR intensities on the number of k points involved in the discretization 
of the Brillouin zone with HF and DFT (see Table IITT1) . As a representative DFT method 
we used the non-empirical generalized gradient approximation (GGA) functional of Perdew- 
Burke-Ernzerhof 31 (PBE). These model calculations suggest that 100 - 1000 k points are 
enough to obtain IR intensities converged within 1%. A large difference in sensitivity of 
vibrational frequencies and IR intensities to should not be considered as a serious issue, 
since the overhead from enlarging Nk is negligible in our vibrational frequency evaluation 
algorithm. 4 We suppose that the large difference in magnitudes of HF and PBE IR intensities 
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TABLE II: Comparison of dipole derivatives with respect to 
atom (in a.u.) for periodic and oligomeric calculations of the 



Cartesian coordinates of the oxygen 
water chain with HF (see Fig. [1]) . 



Method 



8d z /dO z 



dd z /do y 



ddy/dO z 



ddy/dOy 



Oligomeric n/n — 1 

32/31 -1.02563 

52/51 -1.02580 

82/81 -1.02585 

102/101 -1.02587 

Numerical (PBC) 

N k = 1 000 -1.02611 

N k = 10 000 -1.02611 

Analytic (PBC) 

N k = 1 000 -1.02581 

N k = 10 000 -1.02590 

Numerical (PBC) 

N k = l 000 -0.83322 

N k = 10 000 -0.83322 
Analytic (PBC) 

N k = 1 000 -0.79599 

N k = 10 000 -0.79627 



Basis 6-3 1G 

-0.08673 
-0.08675 
-0.08676 
-0.08676 

-0.08664 
-0.08664 

-0.08670 
-0.08676 
Basis 6-311++G** 

-0.09828 
-0.09826 

-0.10754 
-0.10774 



-0.12885 
-0.12885 
-0.12885 
-0.12885 

-0.12872 
-0.12872 

-0.12885 
-0.12885 



-0.08772 
-0.08772 

-0.08496 
-0.08496 



-0.51942 
-0.51942 
-0.51942 
-0.51942 

-0.51931 
-0.51931 

-0.51942 
-0.51942 



-0.54018 
-0.54018 

-0.53669 
-0.53669 



is related to the tendency of pure DFT functionals to yield a more metallic description than 
the HF method. 



B. Poly(paraphenylenevinylene) 

In order to illustrate the prediction capabilities of our approach, we chose the experimen- 
tally and theoretically well-studied one-dimensional system of poly(paraphenylenevinylene) 
(PPV, Fig. [2]) J & i & M i M Such an interest is motivated by conductivity and luminescent nonlin- 
ear properties which make this system a promising candidate as a material for light-emitting 



11 



TABLE III: Convergence of vibrational frequencies (cm -1 ) and IR intensities (km/mol) with the 
number of k points (Nk) for the five normal modes of the water chain in HF and PBE (see Fig. [T|). 
Band gaps are denoted by E g . 





1 






2 




3 


4 




5 




N k 


Freq. 


Int. 


Freq. 


Int. 


Freq. 


Int. 


Freq. 


Int. 


Freq. 


Int. 










HF/6-311++G**, 


Eg = 14.4 eV 








30 


399.0 


92.0 


585.8 


271.7 


1785.3 


146.2 


4069.0 


208.5 


4211.9 


79.2 


50 


399.0 


92.8 


585.8 


271.7 


1785.3 


146.6 


4069.0 


210.8 


4211.9 


79.7 


100 


399.0 


93.3 


585.8 


271.7 


1785.3 


147.0 


4069.0 


212.6 


4211.9 


80.1 


1000 


399.0 


93.7 


585.8 


271.7 


1785.3 


147.4 


4069.0 


214.1 


4211.9 


80.4 


10 000 


399.0 


93.8 


585.8 


271.7 


1785.3 


147.5 


4069.0 


214.2 


4211.9 


80.5 










PBE/6-311++G**, E a = 5.3 


eV 








30 


503.6 


75.8 


544.6 


192.7 


1630.9 


157.9 


3489.0 


426.6 


3789.1 


34.8 


50 


503.6 


76.7 


544.6 


192.7 


1630.9 


157.2 


3489.0 


446.5 


3789.1 


35.0 


100 


503.6 


77.3 


544.6 


192.7 


1630.9 


156.9 


3489.0 


461.4 


3789.1 


35.0 


1000 


503.6 


77.8 


544.6 


192.7 


1630.9 


156.6 


3489.0 


474.6 


3789.1 


35.0 


10 000 


503.6 


77.9 


544.6 


192.7 


1630.9 


156.6 


3489.0 


475.9 


3789.1 


35.0 



diodes. 32 We would like to point out that previous theoretical work from our group on PPV 
was done with some limitations which are overcome in the current study. In Ref. llCj, the IR 
intensities took into account only the changes in the non-periodic part of the dipole moment 
(d\). This caused an underestimation of IR intensities for those vibrations that modify the 
longitudinal dipole moment of PPV. 

As mentioned in the Introduction, IR intensities usually receive less attention and their 
values are presented only graphically in arbitrary units.— ^ Therefore, in Table HVl we com- 
pare IR-active calculated vibrational frequencies with their experimental counterparts. Note 
that in the case of IR intensities only a qualitative comparison with experimental data can be 
done (see Fig. [3]). In our calculations we used the following functionals: the local spin den- 
sity approximation 3 ^ (LSDA), the meta-GGA of Tao-Perdew-Staroverov-Scuseria^ 6 - (TPSS), 
TPSS hybrid 36 (TPSSh), and long-range corrected hybrid PBE^I (LC-a>PBE). Our choice 
was motivated by a previous study of the performance of the TPSS functional in molecules 
and solids.-! 3 ^ We also included the LC-cuPBE functional because we expect some overesti- 
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mation of IR intensities by regular functionals due to the well-known problem with electric 
field response properties in extended systems.—^ According to Table HVl and Fig. [33 all func- 
tionals perform quite adequately, although they generally underestimate lower frequencies 
and overestimate higher ones. Mean absolute errors with respect to experimentally observ- 
able frequencies indicate that the pure DFT methods (LSDA and TPSS) describe vibrational 
frequencies better than the hybrid functionals (TPSSh and LC-cjPBE). Calculated values of 
IR intensities qualitatively follow the right trends in most of the cases. It is remarkable that 
even though LC-cjPBE does not perform very well for vibrational frequencies, it predicts IR 
intensities which are in a better agreement with experiment than are those from the regular 
functionals. This complies with the idea that IR intensities with regular functionals suffer 
mostly from the wrong asymptotic behavior of the exchange potential. It is possible that 
the correct 1/r asymptotic behavior of the LC-cuPBE exchange potential^! is responsible for 
the improvements of IR intensities in elongated systems. 



IV. FINAL REMARKS 

We have presented a simple route to evaluating IR intensities in solids within the HF 
and DFT frameworks with localized basis sets. As in the molecular case, to evaluate IR 
intensities one needs to differentiate the dipole moment with respect to nuclear coordinate 
displacements. However, the periodic dipole moment is not a straightforward generalization 
of its molecular counterpart but rather can be seen as a geometric quantum phase. Thus, we 
have used the discretized Berry phase expression for the dipole moment per unit cell devel- 
oped in the modern theory of polarization^ and elaborated for localized basis sets by Kudin 
and coworkers.— We have differentiated the discretized dipole moment expression with re- 
spect to in-phase nuclei displacements and have demonstrated validity of our technique on 
the model one- dimensional chain of water molecules. The evaluation of IR intensities in- 
troduces only a negligible overhead in the characterization of vibrational frequencies, since 
CPU timings for dipole derivatives [Eqs. (1251) . (1521) . and (1331) ]. in all studied cases constitute 
less than 4 % of total CPU times. The PPV study with different functionals reveals that 
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TABLE IV: PPV IR active vibrational frequencies (cm 1 ) and their intensities (km/mol) calculated 
with various methods and the 6-31G** basis set. In all calculations Nk = 1000 was used. 



Mode 6 




LSDA 




TPSS 


TPSSh 


LC-cjPBE 


Expt. a 


Freq. 


Int. 


Freq. 


Int. 


Freq. 


Int. 


Freq. 


Int. 


Freq. 


A u 


224 


0.2 


225 


0.2 


230 


0.2 


241 


0.2 




B u 


422 


26.5 


421 


23.6 


426 


19.8 


431 


7.1 


429 


A-u 


552 


11.3 


555 


9.7 


564 


10.4 


584 


12.2 


555(s) c 


B u 


795 


30.8 


788 


31.0 


797 


27.2 


809 


12.4 


785 




818 


27.9 


830 


31.4 


843 


32.9 


877 


32.1 


837(s) 




911 


0.4 


941 


0.2 


959 


0.2 


1009 


1.0 




A-u 


948 


33.0 


981 


29.9 


997 


32.4 


1023 


37.5 


965(s) 


B u 


1000 


4xl0~ 2 


1011 


3xl0~ 2 


1025 


2xl0~ 2 


1050 


3.2 


1013 


B u 


1098 


8.5 


1122 


6.3 


1137 


5.9 


1159 


7.0 


1108 


B u 


1200 


14.0 


1236 


2.0 


1247 


1.7 


1240 


1.1 


1211 


B u 


1302 


17.3 


1300 


9.5 


1312 


9.9 


1 1 A 

1314 


1 A O 

14.3 


1271 


B u 


1401 


63.5 


1389 


40.2 


1401 


34.4 


1386 


15.5 


1339 


B u 


1467 


1.5 


1454 


7.8 


1470 


8.1 


1497 


8.1 


1423(s) 


B u 


1538 


109.6 


1539 


81.0 


1561 


79.3 


1608 


61.8 


1518(s) 


B u 


3061 


135.4 


3121 


114.3 


3156 


104.3 


3220 


62.7 




B u 


3090 


60.9 


3136 


84.9 


3171 


77.5 


3230 


38.3 




B u 


3109 


20.2 


3162 


49.7 


3197 


45.1 


3256 


21.1 




MAE d 


21 




17 




27 




43 







33 and 



"References 

^The normal modes are classified according to the crystallographic point group Cik- 
c "(s)" stands for strong banks. 

d Mean absolute error with respect to experimentally observed frequencies. 



although the calculated vibrational frequency adequately reproduce those from experiment, 
the corresponding IR intensities do not always follow qualitatively correct trends. Appli- 
cation of the LC-cjPBE functional corrects IR intensities presumably due to the right 1/r 
asymptotic behavior of the LC-c<jPBE exchange potential but worsens frequencies. We hope 
that our scheme and the results reported here will stimulate further functional development 
in the future. 
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FIG. 1: Schematic structure of the one-dimensional water chain. Optimized parameter values are 
listed in Tabled 




FIG. 2: Poly(paraphenylenevinylene) (PPV). 
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FIG. 3: IR spectra of PPV. The experimental spectrum (bottom graph) is based on data given in 
Ref. 34. The calculated spectra were obtained by using a Lorentzian broadening with a 10 cm -1 



width. 
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